#include "field.hpp"
#include "meshTypeFunctions.hpp"
#include "bc_freeOutflow_velocity.hpp"
#include "bc_freeOutflow_waterdepth.hpp"
#include "util.hpp"

#include <catch2/catch_test_macros.hpp>
#include <catch2/catch_approx.hpp>

using namespace hms;
using Type = hms::boundaryCondition::Type;
using Catch::Approx;

TEST_CASE( "FreeOutflow_velocity" ){

	auto test = []( auto&& dummyMesh ){
		using MeshType = remove_qualifiers<decltype(dummyMesh)>;
		auto mesh { uniMeshAs<MeshType>(3, 4, 1, 1, true) };

		using BC = FreeOutflow_velocity<MeshType>;

		const BoundaryPatch& left   { mesh.boundaryPatch("left") };
		const BoundaryPatch& right  { mesh.boundaryPatch("right") };
		const BoundaryPatch& bottom { mesh.boundaryPatch("bottom") };
		const BoundaryPatch& top    { mesh.boundaryPatch("top") };

		Field f { mesh, 2};

		auto setPatchesToBC = [&](){
			for ( const BoundaryPatch& p : mesh.boundaryPatches() )
				f.setBoundaryCondition( BC{f, p} ).update();
		};

		SECTION("Flow from lower left corner to upper right corner") {

			f.values() = 1.0;
			setPatchesToBC();

			REQUIRE( (f.ghostValues(left)   == -1.0).all() );
			REQUIRE( (f.ghostValues(right)  ==  1.0).all() );
			REQUIRE( (f.ghostValues(bottom) == -1.0).all() );
			REQUIRE( (f.ghostValues(top)    ==  1.0).all() );

			REQUIRE(
				f.boundaryCondition(left).type() == Type::freeOutflow_velocity
			);
		}

		SECTION("Flow from upper right corner to lower left corner") {

			f.values() = -1.0;
			setPatchesToBC();

			REQUIRE( (f.ghostValues(left)   == -1.0).all() );
			REQUIRE( (f.ghostValues(right)  ==  1.0).all() );
			REQUIRE( (f.ghostValues(bottom) == -1.0).all() );
			REQUIRE( (f.ghostValues(top)    ==  1.0).all() );
		}
	};
	SECTION("Mesh"){
		test(Mesh{});
	}
	SECTION("UniMesh"){
		test(UniMesh{});
	}
	SECTION("RectMesh"){
		test(RectMesh{});
	}
	SECTION("StructMesh"){
		test(StructMesh{});
	}
}

TEST_CASE( "FreeOutflow_waterdepth" ){
	auto test = []( auto&& dummyMesh ){
		using MeshType = remove_qualifiers<decltype(dummyMesh)>;
		auto mesh { uniMeshAs<MeshType>(3, 4, 1, 1, true) };

		using BC = FreeOutflow_waterdepth<MeshType>;

		const BoundaryPatch& left   { mesh.boundaryPatch("left") };

		Field d { mesh, 1};
		scalar thresholdOne { 5.0 };
		scalar thresholdTwo { 15.0 };
		scalar waterdepthDomain { 10.0 };

		auto setBC = [&]( scalar limit ){
			d.setBoundaryCondition( BC{d, left, true, limit} ).update();
		};
		
		SECTION("Limited. All waterdepths above the threshold"){
			d.values() = waterdepthDomain;

			setBC(thresholdOne);

			REQUIRE( (d.ghostValues(left) == waterdepthDomain).all() );

			REQUIRE(
				d.boundaryCondition(left).type() == Type::freeOutflow_waterdepth
			);
		}
		
		SECTION("Limited. All waterdepths below the threshold"){
			d.values() = waterdepthDomain;

			setBC(thresholdTwo);

			REQUIRE( (d.ghostValues(left) == thresholdTwo).all() );
		}

		SECTION("Limited. Some above, some below threshold"){
			d.values() = waterdepthDomain;
			d.domainValues(left)(0, 1) = 1.0;

			setBC(thresholdOne);

			REQUIRE( d.ghostValues(left)(0, 0) == waterdepthDomain );
			REQUIRE( d.ghostValues(left)(0, 1) == thresholdOne );
			REQUIRE( d.ghostValues(left)(0, 2) == waterdepthDomain );
			REQUIRE( d.ghostValues(left)(0, 3) == waterdepthDomain );
		}
	};

	SECTION("Mesh"){
		test(Mesh{});
	}
	SECTION("UniMesh"){
		test(UniMesh{});
	}
	SECTION("RectMesh"){
		test(RectMesh{});
	}
	SECTION("StructMesh"){
		test(StructMesh{});
	}
}